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It is generally believed that the angular resolution of the Laser Interferometer Space Antenna 
(LISA) for binary supermassive black holes (SMBH) will not be good enough to identify the host 
galaxy or galaxy cluster. This conclusion, based on using only the dominant harmonic of the binary 
■ SMBH signal, changes substantially when higher signal harmonics are included in assessing the 

' parameter estimation problem. We show that in a subset of the source parameter space the angular 

, resolution increases by more than a factor of 10, thereby making it possible for LISA to identify the 

^S) ■ host galaxy/galaxy cluster. Thus, LISA's observation of certain binary SMBH coalescence events 

could constrain the dark energy equation of state to within a few percent, comparable to the level 
O ' expected from other dark energy missions. 
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I. INTRODUCTION 



2 ' An outstanding issue in present day cosmology is the physical origin of dark energy (see, e.g., Ref. [l| for a review). 
' Probing the equation-of-state-ratio (w(z)) provides an important clue to the question of whether dark energy is truly 
\ a cosmological constant (i.e., w = —1). Assuming the Universe to be spatially flat, a combination of the Wilkinson 
'— 'i Microwave Anisotropy Probe (WMAP) and Supernova Legacy Survey (SNLS) data yields significant constraints on 
' w = — 0.967j;Q Q72 0- Without including the spatial flatness as a prior into the analysis, WMAP, large-scale structure 
^ , and supernova data place a stringent constraint on the dark energy equation of state, w = —1.08 ± 0.12. The Laser 
' Interferometer Space Antenna (LISA) could play an important role in investigating the nature of dark energy as 
fS| ', argued in Refs. [1, 

0^ ' Binary supermassive black holes (SMBH), often referred to as gravitational- wave (GW) "standard sirens" (analogous 
CO I to the electromagnetic "standard candles") 01, are potential sources for the planned LISA mission. LISA would be 
i able to measure the "redshifted" masses of the component black holes and the luminosity distance to the source with 

■ good accuracy for sources up to redshifts of a few. However, GW observations alone cannot provide any information 
C , about the redshift of the source. If the host galaxy or galaxy cluster is known one can disentangle the redshift from the 

' masses by optical measurement of the redshift. This would not only allow one to extract the "physical" masses, but 
^ , also provide an exciting possibility to study the luminosity distance-redshift relation providing a totally independent 
k>( ' conflrmation of the cosmological parameters. Further, this combined observation can be used to map the distribution 
rN \ of black hole masses as a function of redshift 0,0, 111- For this to be possible, LISA should (a) measure the luminosity 

■ distance to the source with a good accuracy and (b) localize the coalescence event on the sky with good angular 
■ " " ' resolution so that the host galaxy/galaxy cluster can be uniquely identifled. 

Refs. [1, 3] identified two potential problems in using binary SMBH as standard sirens. Firstly, they found that 
LISA'S angular resolution might not be good enough to identify the source galaxy or galaxy cluster, and that other 
forms of identification would be needed, and secondly, they pointed out that weak lensing effects would corrupt the 
distance estimation to the same level as LISA's systematic error on the measurement of the luminosity distance. 
Their analyses, like most other in the literature on LISA parameter estimation, were based on the so-called restricted 
post-Newtonian (PN) waveforms. The restricted waveforms (RWF) retain only the leading order (i.e., Newtonian) 
term in the wave amplitude, a PN series, but incorporate the phase up to the maximum available PN order, which 
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is currently 3.5PN [l^, [HI, Recent studies have shown that the inclusion of higher order amplitude terms in 
the waveform (and hence higher harmonics of the orbital frequency) would play an important role in the detection 
rates (by increasing the mass reach of the detector) [l^, [3, EH as well as in the problem of parameter estimation 
[l6l. [TtI. Iisl. [l9l. [20j of both ground-based and space-based detectors. Specifically, Refs. [3, [3 examined the improved 
angular resolution of different space-based detector configurations due to the inclusion of higher harmonics. 

In the present work, we revisit the problem of parameter estimation in the context of LISA using amplitude- 
corrected PN waveforms. We investigate systematically the variation in parameter estimation with PN orders by 
critically examining the role of higher harmonics in the fast GW phasing, higher PN corrections in the amplitudes 
and frequency sweep and their interplay with the slow modulations induced due to LISA's motion. More importantly, 
we explore the improvement in the estimation of the luminosity distance and the angular parameters due to the 
inclusion of higher harmonics in the waveform. We translate the error in the angular resolution to obtain the number 
of galaxies (or galaxy clusters) within the error box on the sky. We find that independent of the angular position of 
the source on the sky, higher harmonics improve LISA's performance on both counts raised in Refs. 0, On the one 
hand we will show that the angular resolution improves typically by a factor of ~ 2-500 (greater at higher masses) 
and the error on the estimation of the luminosity distance goes down by a factor of 2-100 (again, larger at higher 
masses). For many possible sky positions and orientations of the source, the inaccuracy in our measurement of the 
dark energy would be at the level of a few percent, so that it would only be limited by weak lensing. We conclude that 
LISA could provide interesting constraints on cosmological parameters, especially the dark energy equation-of-state, 
and yet circumvent all the lower rungs of the cosmic distance ladder. 

This paper is structured as follows. In the next section we introduce our signal model and the LISA noise power 
spectral density we will use. In Section III we discuss our results on parameter estimation and their relevance for 
astrophysics and cosmology. Section IV gives an overview of various effects that are likely to affect our estimates. 
Conclusions arc presented in Section V. Technical details on how the parameter estimation was performed can be 
found in Appendix A. Finally, in Appendix B we give an in-depth discussion of the way parameter estimation is 
influenced by the inclusion of higher harmonics and their amplitude corrections. 



II. SIGNAL MODEL AND LISA NOISE POWER SPECTRAL DENSITY 



The post-Newtonian formalism has been used to study the evolution of a binary under gravitational radiation 
reaction to a very high order in the small parameter v characterizing the velocity of the component objects, yielding 
accurate expressions for the orbital phase and the two gravitational wave polarizations. For binaries consisting of 
component stars of negligible spin on quasicircular orbits, the most accurate computations currently known have 
corrections not only to the orbital phase up to order v"^ (i.e., 3.5PN order in the notation of PN theory^ 
HH, H^, but also corrections to the gravitational wave polarizations to order (i.e., 2.5PN order) (23l. \24 . 
shall call this the "full" waveform (FWF). 

The waveform as seen in LISA is modulated in two ways due to LISA's motion. LISA consists of three spacecraft 
at the vertices of an equilateral triangle of 5 million kilometers, each craft on a heliocentric orbit slightly inclined 
to the ecliptic. As the craft orbit the Sun, the triangular formation also spins around itself with the same one-year 
period as the orbital period. Therefore, relative to LISA the source location and orientation changes with time with 
a one-year period and induces amplitude and phase modulations in the waveform. 

It is well-known that at signal frequencies / 5 x 10""^ Hz, LISA can essentially be modeled as a pair of two-arm 
interferometers, usually labeled as I and II |26| . and this suffices for the sources considered in this paper. (However, 
it would be interesting to investigate the added value of the remaining third combination ignored in this work.) In 
what follows, to begin with we consider a single detector. 

Let us consider a source of total mass M = mi -\- TO2 and symmetric mass ratio v = 17111712/ M'^ (where mi, TO2 
are the individual component masses) located at a luminosity distance Dl- In the stationary phase approximation 
(SPA), the Fourier transform hi( f) of the response of detector I to the full waveform, including the modulations due 
to LISA'S motion, is given by |15j : 

^ ^ A{,^„^,)(t(/.)).t+i(t(A.))e-'^;^w.)W/0) 
hiif) = -:r~n~ 1 exp [i 7/./,^, (*(/&))] , (2.1) 

k=ln=0 2^kF{t{fk)) 

where fk = f /k, an overdot denotes derivative with respect to time, and ipf.kitifk)) is given by 

^/.fc(t(A.)) = 27r/i(A) -k^{t{fk)) - kMtifk)) - 7T/4. (2.2) 

In Eq. (j2.ip . the factor ■\/3/2 is due to the 60° angle between the interferometer's arms and t{f) is given in Ref. [ll|. 
The waveform is a superposition of harmonics of the orbital frequency (labeled by the index k), and each harmonic has 
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PN contributions to the amplitude (labeled by n; note that we can only go up to n = 5, as no amplitude corrections 
are explicitly known beyond 2.5PN). As the PN order in amplitude is increased, more and more harmonics appear; at 
2.5PN order there are seven, which is why the index k only runs up to = 7. Quantities in Eqs. p.ip and (|2.2p with 
the argument t{fk) denote their values at the time when the instantaneous orbital frequency sweeps past the value 
f /k and x{t) is the PN parameter given by x{t) = (27rA/F(t))^/'^, F{t) being the instantaneous orbital frequency of 
the binary. Aj^, „y2)(^) ^^^^ "^^Ifc n/2)(^) '"^^^ polarization amplitudes and phases of the fcth harmonic appearing at 
the n/2th PN order, ^'(t) is the orbital phase of the binary and (?!)d(0 is a time-dependent term representing Doppler 



modulation. Explicit expressions for ^|^. „y2) ^^"^ '^\k n/2) '-■^^ ^'^ found in time-dependence of these quantities 
arises through the beam-pattern functions due to the varying sky position and orientation of the source relative to 
the detector The expression for (t>Dit) is given in [l^,l2^. For the PN expansions for t{F), '^{F), F{F) we refer 
to Ref. 

The restricted waveform (RWF) corresponds to retaining the term with fc = 2 and n = in Eq. (|2.ip and neglecting 
all others. It is clear that the RWF has only the dominant harmonic at twice the orbital frequency but no other 
harmonic, nor PN corrections to the dominant one. It does, however, include the post-Newtonian expansion of the 
phase to all known orders, i.e., up to w^. The FWF, on the other hand, has not only the dominant harmonic but 
also other harmonics up to seven times the orbital frequency and their PN corrections to the relevant order. The 
distinctive nature of the FWF as compared to RWF, especially the richer structure in its spectrum, can be seen in 
Fig. 1 of Refs. 0,[23l. 

Expanding the denominator and extracting the lowest order term helps us rewrite the waveform in a form more 
suitable for our purposes, 

8 7r2/3 i2fy/'t[ 



\n=0 ) \m=0 



(2.3) 



where M. — M v^l^ is the chirp mass of the binary, and [ • ]p denotes consistent truncation to pth post-Newtonian 
order (in our case p = 2.5). The coefficients S'(m/2) are the PN expansion coefficients of F{t{fk))~^/^ and are given 
in Eq. (A.4) of [H. 

Each harmonic in hi{f) is taken to be zero outside a certain frequency range. The upper cutoff frequencies are 
dictated by the last stable orbit (LSO), beyond which the PN approximation breaks down. For simplicity we assume 
that this occurs when the orbital frequency F{t) reaches Flso = 1/{6'^^^2ttM) - the orbital frequency at LSO of a 
test particle in Schwarzschild geometry in c = G = 1 units^. Consequently, in the frequency domain, the contribution 
to h{f) from the fcth harmonic is set to zero for frequencies above fc-FLso- In determining the lower cutoff frequencies 
we assume that the source is observed for at most one year, and the fcth harmonic is truncated below a frequency 
fci^in, where Fin is the value of the orbital frequency one year before LSO is reached (T5| : 

Fin = ^^(iLSO - Atobs) = ^ ^. (2.4) 

For simplicity the quadrupole formula was used. In the above, <lso and wlso = l/%/6 are, respectively, the time 
and orbital velocity at last stable orbit, and Aiobs = 1 yr- However, LISA's sensitivity becomes poorer and poorer 
below ~3 mHz and current estimates normally assume a "noise wall" at /g = 10^"' Hz. Thus, we take the lower cutoff 
frequency of the fcth harmonic to be the maximum of /g and fci^in- For more details we refer to p^ . 

As we mentioned before, the LISA detector can be viewed as a combination of two independent detectors. Until now 
we have dealt with one detector. Calculations for the two detector case closely follow the corresponding treatment 
for the RWF, which is explained in detail in [2^1. The beam-pattern functions for detector II can be obtained 
from those of detector I by a simple rotation. The waveform hu{f) for the second detector is formally identical to 
hi{f), and quantities ^^^-^ and ^^2) ^^''^ obtained from their counterparts ^^2) ^^^^ 'f'lk n/2) replacing the 
beam-pattern functions of detector I by those of detector II. 



^ Note that the cutoffs arc placed on the orbital frequency of the binary, not the dominant GW harmonic; hence the extra factor of 2 in 
the denominator of the expression for -Flso- 
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The waveform given in Eq. (|2.3p together with its counterpart for the second detector form the basis of the analysis 
in this paper. FoUowing earUer works of Refs. [1, 0, [H, [13, [S] we employ the Fisher matrix approach [2^, [s^l to the 
problem of parameter estimation. The waveforms depend on nine parameters which are chosen to be 

p = (lnA^,(S,<c,'^c,ln£'L,Ms,ML,(?!'s,0L) , (2.5) 

where S = (m2 — nii)lM (mi and m2 being the individual masses; throughout this paper we assume m2 > is used 
as a parameter instead of v following Ref. [l^l; ic, 0c are, respectively, the time and orbital phase at coalescence^; 
/is = cos 0s and <j)s determine the source position in the sky; while /iL = cos 0l and (pi^ determine the orientation 
of the binary's orbit with respect to a nonrotating detector at the solar system barycenter"'. Following Ref. [1^, we 
have fixed the initial position and orientation of LISA by setting the constants (/)o and ao defined there to zero at 
i = 0. The polarization amplitudes and phases depend on the location and orientation of the source through the 
beam-pattern functions. The orientation of the source changes relative to LISA with the period of a year. Therefore, 
these quantities are modulated on a one-year time scale and depend on /^s. Ml, 0s, 0l and also on 5, InA^ and Iq. 
The phase ipf,k{t{fk)) depends on In A^, (5, Iq and 0c and varies with the orbital timescale which is much shorter than 
a year. 

Following Ref. [13], the parameters can be subdivided into two subcategories: pfast ~ (In Al, i5, tc, 0c) and Psiow = 
(IuDl, Ms, Ml, 0s, 0l)- The classification arises naturally because the signal that LISA observes can be viewed as a 
(slow) low frequency modulation due to its motion around the Sun superposed on the (fast) high (GW) frequency 
carrier signal. The accuracy of estimation of Pfast follows from the GW phasing of the binary while that of Psiow from 
the modulations associated with LISA's orbital motion. 

In our analysis, we take the noise power spectral density (PSD) to be that given by Eqns. (2.28)-(2.32) of Ref. 
As usual, the Fisher matrix T for LISA as a whole is simply T = Fj + Fn, where Fiji arc the Fisher matrices 
computed from the waveforms hi ii{f). The parameters used will be the ones listed in Eq. (|2.5p . so that F is a 9 x 9 
matrix. However, the errors in the estimation of /^s and 0s obtained in this way will be converted to a solid angle 
Afis centered around the actual source direction. Following the notation of [28| . 

Ar^s = 2^V(Ams A0s)2- (,5/is '50s)2, (2.6) 

where the second term is the covariance between /^s and 0s. A similar quantity Af^L is used to quantify the error 
estimate in the orientation of the binary's orbit. 

In what follows, whenever it is necessary to consider a specific cosmological model we will assume a flat Universe 
with Hubble constant Hq = 75 kms^^Mpc""'^, matter density Jim = 0.27, and dark energy density fid ~ 0.73, with 

f^Total = f^m + = 1- 

III. THE RESULTS AND THEIR ASTROPHYSICAL IMPLICATIONS 

In this Section we will discuss the extent to which LISA will be able to constrain cosmological parameters by 
observing a binary SMBH with a large signal-to-noise ratio and measuring their parameters, most importantly their 
location on the sky and the luminosity distance. As we shall see, the number of clusters in LISA's angular error box 
reduces dramatically when using the FWF as compared to RWF, thereby enabling us to identify the host galaxy, or 
galaxy cluster. Consequently, it should be possible to measure the dark energy equation of state by combining LISA 
observations with electromagnetic observations. 

We start by outlining the generic features of parameter estimation with higher harmonics; after that we will focus 
on angular resolution and the estimation of luminosity distance, and what these can tell us about cosmology. 



^ The waveform is invariant under interchange of mass labels provided that at the same time, the phasing is shifted by tt; since we use a 
single phasing formula we need an ordering in the definition of <5. The parameter 5 was originally introduced in place of u to ensure that 
the FWF Fisher matrix remains nonsingular on the surface mi = m2. The situation is reversed in the case of RWF: there the Fisher 
matrix becomes singular in the limit <5 — ► whereas it remains regular at mi = m2 when u is used in place of 5 (see the discussion in 
[20l | ) . The equal mass case is dealt with in more detail in |3l| . 

^ Below we will consistently set the values of tQ and ipc to zero, but both parameters are included as coordinates on the space of signals 
in computing the Fisher matrix. 

* This is a different notation from Ref. [Ssl , where the source angles measured in the fixed barycenter frame arc denoted by (6s , 0s , , 0L ) ■ 
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A. Parameter estimation with the full waveform: General trends 

Inclusion of higher harmonics results in a significant improvement in the determination of a binary's parameters 
in the context of LISA, as is the case for ground-based detectors [2^. A typical variation of parameter estimation 
(PE) accuracy with PN orders in amplitude is displayed^ in Table [J We observe the following general features 
of amplitude-corrected waveforms with regard to PE: For all masses and all angles we have explored, there is a 
significant improvement in the estimation of all parameters for LISA when considering the full waveform as compared 
to the restricted PN waveform. The orbital frequency and the inspiral rate, and therefore the phase evolution of 
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TABLE I: Variation of parameter estimation errors with post-Newtonian orders in amplitude for a (10®, 1O''')M0 binary at 
z — 0.55 (corresponding to a luminosity distance of Dl ~ 3 Gpc for a Hubble constant Ho ~ 75kms-^Mpc, matter density 
Qm = 0.27, and dark energy density Qd ~ 0.73, with S^Totai ~ + ~ 1). The angles are chosen to be /ig — —0.8, 0s = 1, 
fiL = 0.5, 0L = 3. 

the waveform, are determined primarily by the chirp mass A4. Thus accurate phase tracking leads to a precise 
measurement of M. The phase also depends on the "mass difference" 6 which, therefore, can also be measured quite 
accurately. This is borne out by Table [J where the trend is shown to be true for the FWF also. Finally, note the 
spectacular improvement in angular resolution and the determination of luminosity distance in going from RWF to 
FWF, which will be the focus of the rest of this Section. 

We relegate to Appendix B a more critical and in-depth discussion of the trends observed with increasing PN order 
in amplitude and consequent inclusion of higher harmonics in the waveform. 

B. Effect of higher harmonics on angular resolution and luminosity distance 

Table |ll] lists the one-sigma errors in parameters of interest for seven different combinations of the angular param- 
eters, as in Ref. [l^, each for three different binary masses. The angular parameters (cf. Al, . . . , AT in the Table) 
chosen are a coarse sample of the possible orientations of the orbit and the source's sky location and our choice of 
masses is indicative of the different binary SMBH coalescences LISA is likely to observe with a large signal-to-noise 
ratio (SNR). To make direct contact with astrophysical systems, we give the physical masses M = Mphys and not the 
observed (i.e., redshifted) masses Mobs- The two are related by Mobs = (1 + ^)Mphys, where z is the cosmological 
redshift of the source. In order to compute the upper frequency cutoff one should first convert the total physical 
mass to total observed mass and then use the formula for the LSO frequency. Our sources are all at z ~ 0.55, i.e., 
a luminosity distance of Di^ = 3 Gpc. Theoretical predictions of event rates for SMBH coalescence vary over a wide 
range, but the rate could be as high as 1 per year within z = 0.55 [11]. (See Section IV for a discussion of our results 
for merger events occurring at a higher redshift.) The physical masses and the corresponding LSO frequencies in the 
form (mi/Me, mz/Mo, 7FLso/mHz), are (lO^, 10^ 9.03), (6.45 x 10^ 1.29 x 10^ 7.33) and (10^ 10^ 0.903). Thus, 
the highest harmonic (at seven times the orbital frequency) of the heaviest system that we consider does not quite 
reach the sweet spot of EISA's sensitivity. For the other two systems the dominant harmonic is close to the sweet 
spot and higher harmonics sweep through EISA's sensitivity bandwidth. 

The table lists the 1-sigma errors incurred in the estimation of all the parameters except for the errors on /iL, 0l 
(the direction of orbital angular momentum) and 0c- As usual, we have converted the error in the estimation of 



^ The numerical values in all our tables and results are unaffected by the missing terms discussed in |25|| to the accuracies quoted. 
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and (ps to an error in the solid angle Af^s centered around the actual source direction. For the sake of completeness 
we have given, in alternate rows, the errors for both the RWF and FWF. 

For certain values of the angular parameters the presence of the harmonics seems to have a considerable impact 
on the determination of the luminosity distance and the angular position of the source. The errors in the luminosity 
distance (i.e., ADl) and the source's sky position (i.e.. Ails) are reduced by factors up to 600 and 400, respectively, 
while using FWF as compared to RWF. This means that the error box could be smaller by a factor of 2.4 x 10^. 
However, what is relevant for cosmological applications is by how much the error in the sky position goes down (i.e., 
about a factor of 2.5 to 400). Interestingly, the heavier, and astrophysically most relevant, system of (10®, 10^) Mq, 
where none of the harmonics get close to the detector's sweet spot, shows the largest improvement in distance 
estimation and angular resolution in going from restricted to amplitude-corrected waveforms. This observation, based 
only on the few systems studied in this paper, is found to be generally true in an independent and more exhaustive 
study by Trias and Sintes (3l| . The larger improvement over the RWF is to be expected since for very massive systems 
only the higher harmonics radiate significantly within the detector's bandwidth. However, what is striking is that 
when considering only FWF, for most choices of angles the errors on distance and the angular resolution are almost 
at the same level as for lighter systems. Doppler modulation does not seem to affect the accuracy of estimation of 
parameters for the systems considered in this paper; they are expected to be important for systems with lower masses 
[26|. An alternative method to study these issues is by a direct use of the three time delay interferometry variables. 
Work along these lines is in progress [s^ l and should provide an independent check of our results in the near future. 



C. Number of clusters in LISA's error box 

Binary black holes are standard sirens Q . The amplitude of gravitational waves from a binary SMBH is proportional 
to A4^^^/Di^. As evidenced by the numbers in Table HIl LISA will measure both the chirp mass A4 of the source and 
the amplitude of gravitational waves to a great precision. Thus, the luminosity distance to a source can be extracted by 
gravitational-wave observations alone. In order to derive the luminosity distance-redshift relation, it is also necessary 
to measure the redshift z to the source, but LISA cannot measure z. However, it might be possible to determine 
the source's redshift if the host galaxy, or galaxy cluster, can be optically identified. Whether or not this is possible 
depends on how good LISA's angular resolution is, and whether it is small enough that no more than a few galaxies 
or galaxy clusters are found within the angular error box®. 

To determine the number of galaxy clusters within a solid angle AfJg of the angular error box we need the comoving 
volume corresponding to a cone defined by Afis whose height is the physical distance from the detector to the source, 
which, of course, depends on the cosmological model. In a Universe whose matter density is fim and in which dark 
energy density takes the form of a cosmological constant'' f2d (with f^xotai = + ~ 1), the comoving volume per 
unit redshift within a box of angular size Af2s is: 

dVcjz) _ Arts Dl{z) 1 

dz Ho (1 + zf y/n^^ii + zf + ' 

In the above, Di^{z), the luminosity distance as a hmction of redshift, is given by 

, , 1 + z r dz' , ^ 

Dh(z) = —rr- / , , (3.2) 

^ ' Ho Jo Vf^m(l + z'f + r!d ^ ' 

where Ho is the Hubble parameter at the current epoch. The comoving volume from the observer to the source within 
a cone defined by AOg is simply the integral of Eq. (|3.ip : 



^, Ans Dl{z') 1 

Ho (1 + z'Y va„(i + z'f + 



Vciz) = I dz' -JL^-JA^^===== . (3.3) 



The number density of clusters at high redshifts is not known very well. Following Ref. [3a | we assume that the 
number density of clusters is ~ 2 x 10~^/i^Mpc~'^, where h is the present value of the Hubble parameter in units 
of lOOkms^^Mpc^^. We take h ~ 0.75, 51,„ — 0.27 and fid = 0.73. For this choice of cosmological parameters the 



Note that we have to use the angular error box and cannot use the smaller volume error box also fixed by the luminosity distance. In 
order to do precision cosmology we have to measure the source's redshift independently of the luminosity distance. 
In the next subsection we will relax this assumption. 
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FWF 


2497 


0.17 


26 


2.7 


9.7 


1.1 


0.53 


0.0096 


A4 


-0.5 


3 


-0.6 


-2 


RWF 


1868 


0.74 


150 


3.1 


15 


1.2 


3.1 


- 












FWF 


1781 


0.19 


13 


2.5 


12 


0.58 


0.27 


0.011 


A5 


0.9 


2 


-0.8 


5 


RWF 


3740 


15 


84 


2.3 


8.0 


2.1 


1.7 


0.82 












FWF 


2857 


0.11 


8.1 


1.7 


7.9 


0.69 


0.17 


0.0062 


A6 


-0.6 


1 


0.2 


3 


RWF 


2185 


0.42 


220 


3.9 


15 


2.9 


4.5 


- 












FWF 


2108 


0.24 


65 


3.0 


11 


1.6 


1.3 


0.014 


A7 


-0.1 


3 


-0.9 


6 


RWF 


2213 


0.58 


410 


3.5 


13 


1.1 


8.4 


- 












FWF 


2175 


0.45 


300 


2.9 


10 


0.74 


6.1 


- 


(mi, m2) = 


(6.45 X 10'', 1.29 X lO'')i\f0 
















Al 


0.3 


5 


0.8 


2 


RWF 


385 


1.3 


21 


5.5 


13 


3.2 


0.43 


0.073 












FWF 


511 


1.0 


8.4 


4.2 


9.1 


2.1 


0.17 


0.056 


A2 


-0.1 


2 


-0.2 


4 


RWF 


595 


1.1 


120 


4.2 


9.2 


2.5 


2.4 


0.062 












FWF 


771 


0.70 


25 


3.3 


6.5 


1.7 


0.51 


0.039 


A3 


-0.8 


1 


0.5 


3 


RWF 


1345 


0.33 


170 


3.4 


5.8 


2.7 


3.5 


- 












FWF 


1573 


0.25 


53 


2.6 


4.2 


1.6 


1.1 


0.014 


A4 


-0.5 


3 


-0.6 


-2 


RWF 


924 


0.78 


160 


3.0 


6.8 


1.7 


3.3 


- 












FWF 


1158 


0.26 


27 


2.3 


5.0 


1.0 


0.55 


0.015 


A5 


0.9 


2 


-0.8 


5 


RWF 


1863 


15 


87 


2.4 


3.8 


2.2 


1.8 


1.0 












FWF 


1506 


0.19 


25 


2.0 


3.9 


1.3 


0.51 


0.011 


A6 


-0.6 


1 


0.2 


3 


RWF 


1069 


0.47 


240 


4.1 


7.2 


3.1 


4.9 


- 












FWF 


1378 


0.32 


110 


2.9 


4.8 


2.1 


2.2 


0.018 


A7 


-0.1 


3 


-0.9 


6 


RWF 


1093 


0.57 


420 


3.1 


6.1 


1.6 


8.6 


- 












FWF 


1448 


0.50 


350 


2.5 


4.2 


1.1 


7.1 


- 


(mi, m2) = 


(lo^ 


10'')A/o 






















Al 


0.3 


5 


0.8 


2 


RWF 


495 


11 


600 


1400 


1100 


290 


12 


- 












FWF 


444 


2.2 


16 


190 


240 


75 


0.33 


0.12 


A2 


-0.1 


2 


-0.2 


4 


RWF 


773 


10 


6500 


870 


710 


190 


130 


- 












FWF 


685 


1.2 


43 


130 


160 


51 


0.88 


0.068 


A3 


-0.8 


1 


0.5 


3 


RWF 


1824 


3.2 


24000 


380 


310 


90 


490 


- 












FWF 


1549 


0.29 


82 


96 


100 


31 


1.7 


0.016 


A4 


-0.5 


3 


-0.6 


-2 


RWF 


1249 


6.9 


2400 


550 


450 


120 


49 


- 












FWF 


1081 


0.34 


40 


110 


130 


40 


0.82 


0.019 


A5 


0.9 


2 


-0.8 


5 


RWF 


2493 


110 


8300 


270 


220 


63 


170 


- 












FWF 


1954 


0.18 


18 


200 


180 


49 


0.37 


0.010 


A6 


-0.6 


1 


0.2 


3 


RWF 


1465 


4.7 


53000 


470 


380 


110 


1100 














FWF 


1273 


0.44 


300 


105 


115 


36 


6.1 




A7 


-0.1 


3 


-0.9 


6 


RWF 


1480 


21 


170000 


520 


390 


98 


3500 














FWF 


1300 


1.3 


3400 


87 


100 


30 


69 





TABLE II: Accuracy in LISA's measurement of the various parameters in Eq. (|2.5|) . for seven different sets of the angular 
parameters and three different combinations of the (physical) masses at a distance of 3 Gpc {z = 0.55). When the number of 
clusters in the error box on the sky is significantly larger than 1, it will not be possible to determine redshift unless the inspiral 
event has a clear optical counterpart; we have chosen not to quote results for Aw in such cases. (Note that the error on w is 
ultimately determined by both LISA's statistical errors and weak lensing errors in the determination of luminosity distance.) 
The figures clearly demonstrate significant improvement in parameter estimation when higher order terms are included. 
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number of expected clusters within the LISA angular error box is given in the second last column of Table |TT1 Clearly, 
in many cases the number of clusters is of order 1 which means that LISA will help us to uniquely identify the host 
galaxy cluster of a binary SMBH merger event within a redshift of z = 0.55. 

There is a caveat with regard to the number of clusters found within the error box that is important to mention 
at this stage. Note that we integrated the comoving volume up to z = 0.55, the location of our source. In reality, we 
would not know the redshift of the source and are not fully justified in integrating only up to this point; in principle 
we should consider all clusters within the error box up to much larger redshifts. However, sources at redshifts much 
larger than z = 0.55 for the same luminosity distance would probably give radically different cosmological parameters. 
Consistency with other observations justifies considering only galaxy clusters that are roughly in the redshift region 
determined by inverting a luminosity distance-redshift relation based on parameter values from other measurements. 
In any case, we have checked that integrating the comoving volume up to z = 1 (which for our chosen value of Z3l 
would already imply a very significant departure from current cosmological models) does not drastically change the 
results of Table [III Indeed, despite the higher limit of integration, the number of clusters in the angular error box 
remains less than 3 for most choices of angles, which is the (arbitrary) cutoff in Nciustcis we have chosen to assess 
whether redshift can be determined. 

There have been suggestions that in order to identify the source associated with a binary SMBH merger we should 
also look for optical/UV counterparts; the improved angular resolution with the correct signal model should help in 
this case too. Optical and other electromagnetic telescopes will need to survey a much smaller area on the sky than 
was thought before and should therefore more easily identify the galaxy cluster in which the merger took place. 

The error in Dl being less than a percent means that we should, in principle, be able to tightly constrain the 
cosmological model a lot better. However, as discussed by several groups, the possible effect of weak gravitational 
Icnsing on parameter estimation, in particular on our ability to measure the luminosity distance (see, e.g., Ref. 0]), 
will limit the extent to which LISA can measure dark energy. These considerations do not alter the main conclusions 
of this paper as our main goal is to show that the field of view in LISA's angular resolution in most cases involves 
only a few sources. 



We conclude by mentioning the implications of our estimates for astrophysics and cosmology. Up to this point 
in our analysis we have not assumed any electromagnetic counterpart to the coalescence events LISA will observe. 
But in reality, many of the LISA observations are likely to have electromagnetic counterparts, either as a precursor 
or as an afterglow [36l . [37| . The implications of an uniquely identified quasar source in coincidence with LISA were 
examined in [38| specifically assuming quasars as a possible electromagnetic counterpart. Even a single event of this 
type would provide us with unprecedented tests of SMBH accretion physics, such as a precise measurement of the 
Eddington ratio [sst . For this to be possible, and to identify an unique electromagnetic event in coincidence with a 
merger event as seen by LISA, the use of higher harmonics would be crucial since for many of the directions in the 
sky it brings down the number of candidate clusters to less than one. 

Another exciting possibility is to use LISA as a cosmological probe. If a unique host is identified electromagnetically 
in coincidence with the LISA observation, the redshift of the host galaxy will be known to very high accuracy. The 
improved estimate of the luminosity distance obtained by using the FWF would play a crucial role in determining 
the cosmological parameters as suggested by Schutz Q . 

Gravitational-wave observation of a single inspiral event coupled with an electromagnetic determination of the 
redshift would allow LISA to strongly constrain the equation of state of dark energy. In a flat Universe, the luminosity 
distance can be written as 



where H(z) is the Hubble parameter. Given a form of matter energy with density parameter and a (constant) 
equation-of-state index w, one has 



D. Constraints on the equation-of-state of dark energy 




(3.4) 



H{z) = Ho an(l + zf + fliil + z)3(i+-) 



nl/2 



(3.5) 



For a given, fixed redshift (and fixed Hq, f^m, and ^d), the error on w is 



1 



Aw = 



dw 



(3.6) 
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Using (|3.4p and (|3.5p and setting fl^ = 0.73, w 



1, and Dl = 3 Gpc. one obtains 



533 Mpc. 



(3.7) 



dw 



With these assumptions and using the values for ADi^/Di^ from our analysis, we find that in all of our examples, the 
FWF always leads to a smaller value for Aw than the RWF whenever a comparison can be made. Indeed, in Table 
im there are many instances where the RWF doesn't allow for a determination of the redshift because the number 
of clusters in the angular error box on the sky is too large, in which cases w cannot be measured. In most of the 
examples we have considered, the FWF does not have this problem. 

The foregoing analysis does not take into account the error in luminosity distance arising as a result of weak lensing 
of gravitational waves by the intervening mass concentrations between the binary SMBH source and LISA. It is 
estimated that the weak lensing will introduce errors in the luminosity distance at the level of about 3-5% for sources 
at z ~ 0.5 0]. This is far greater than systematic error in LISA's measurement of the luminosity distance. Therefore, 
weak lensing will be the limiting factor in LISA's ability to measure the dark energy equation of state unless weak 
lensing can be corrected by properly modeling the weak lenses. 

We now revisit the caveat in the previous subsection regarding the fact that we tend to consider an error cone that 
stretches only to z = 0.55. Let us evaluate Aw when considering a measured luminosity distance of Z^l = 3 Gpc, but 
allowing for the possibility that the source may be at a different redshift. Given our luminosity distance, to assume 
that the source could be at, e.g., z = 0.6 and keeping Hq, and Qd the same would already require w = —0.47 
(from Eqs. ([331) and ^M): a value excluded by WMAP and SNLS results 0. In that case \dDL/dw\ ~ 668 Mpc; 
substituting this into the right hand side of Eq. (|3.6p one finds that all values of Aw will decrease by about 25 percent 
compared with the z ~ 0.55 (w = —1) case. In reality, uncertainties in Hq, Jim, and fij will also have to be taken 
into account, but it will probably not be necessary to consider potential sources in the angular error box at redshifts 
that differ from the "favored" value by more than 20 percent. A more complete treatment should take into account 
the covariances among Hq, fly^, ^d, and w] this we relegate to a future study. 

We conclude by noting another interesting feature of our analysis. In addition to the improved angular resolution 
and luminosity distance, the errors in the estimation of mass parameters reduce considerably while using FWF. 
Even though the RWF itself would give a very good estimate of the masses, the improved measurement of the 
mass parameters would be very important in performing certain novel tests of general theory relativity and its 
alternatives 0, [H, H^. For example, the use of FWF should help improve the accuracy with which the individual 
phasing coefficients can be determined independent of each other, an idea proposed in Ref. [39l. l40j. 



In this Section we will briefly discuss some caveats regarding our results. Here we focus on physical issues; limitations 
related to our chosen methods for computing errors (and their resolution) are commented on in Appendix A. 

a. Effect of black hole spins. Our analysis is restricted to the case of nonspinning black holes while astrophysical 
black holes are known to have significant amounts of spin. Including the spin effect could significantly alter the 
estimation of different parameters . To get an estimate of the effect of the spin on parameter estimation, we have 
recomputed the covariance matrix by including the leading order spin-orbit parameter in the waveform model. We 
observed a deterioration of up to an order of magnitude for the mass parameters but only a factor of a few in the 
estimation of luminosity distance and angular resolution. As suggested by Vecchio [4l|, the inclusion of precession 
should compensate for the deterioration of the mass errors to a great extent. For instance, including the leading 
spin-orbit term with a simple precession model, Ref. [4]| showed that the angular resolution and luminosity distance 
determination could be improved by a factor between 3-10. A more recent analysis [i^l incorporated the 2PN spin 
effects. Thus our estimates may not be too far from the realistic case despite the assumption of nonspinning holes. 
An interesting exercise for the future would be to include the amplitude corrections with the spin effects and see the 
extent to which the results of this paper change. 

b. Weak and strong gravitational lensing. An important physical effect we have neglected in our analysis is 
the possibility of gravitational lensing of the signal. Weak lensing by intervening matter distributions distorts the 
gravitational waveform, inducing a systematic error in the estimation of luminosity distance [1]. This could be as 
high as 5-10% for some of the systems, much higher than LISA's systematic error on the luminosity distance. On the 
other hand, strong lensing may improve the angular resolution, as multiple gravitational-wave "images" formed will 
reduce the correlations between various parameters [1^. This improvement could be as high as 100 times for a million 
solar mass (redshifted) binary [31. Thus the net effect of lensing of a high redshift source may not be disastrous. 



IV. FACTORS AFFECTING THE ESTIMATES 
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c. Inaccuracy in the cosmological parameters. In computing the accuracy with which we can measure the dark 
energy parameter w we have assumed that other cosmological parameters, namely Hq^ and f2dj are all known 
precisely and that the redshift z of the host galaxy is measured accurately. In reality, all these quantities will have their 
own statistical errors that must be folded into the analysis and we hope to do so in a future publication. Moreover, it 
might be possible for LISA to measure the luminosity distance-redshift curve by observing many supermassive black 
hole binary coalescences during its mission lifetime and extract all the cosmological parameters. How well LISA might 
be able to do would be a very interesting question that wc would like to reserve for future studies. In any case, the 
fact that higher harmonics enable a greater accuracy in the measurement of the parameters should help. 

d. Event rates and measurement of w. The rate of binary SMBH mergers in the Universe is not known with good 
accuracy. The rate depends on the cosmological model of choice and on structure formation and growth. Current 
models predict rates that vary over two orders of magnitude. For one of the models the rate of binary SMBH merger 
is as large as ~ 10yr~^ within z = 1 or about ~ lyr~^ up to z = 0.55 (see [s^ for an overview). A single source at 
the right point in parameter space at z = 0.55 will be good enough to measure w to a few percent. However, events 
at z = I will be more frequent and consequently a few of them are likely to end up in the favored region of parameter 
space. We have verified that at z = 1, for all choices of masses considered and angles A1-A5 (except for A3 in the 
case of the most massive system), the number of clusters within the error box is less than 10 when using the FWF. 
It should therefore be possible to either directly identify the merger source (by demanding consistency of the inferred 
cosmological parameters with those from other observations), or to locate it by observing an afterglow. Values of Aw 
will be about 70 percent higher than before; however, with repeated observations one would be able to perform a 
statistical analysis along the lines of Thus, merger events at z = 1 also allow for a measurement of w, although 
in this case the microlensing would probably restrict our accuracy in its estimation to ^ 10 percent 0]. 

V. CONCLUSIONS 

In this paper we have reassessed LISA's ability to perform precision cosmology by using waveform models with far 
richer structure than those used before, and we find remarkable improvement in its angular resolution and estimation 
of the luminosity distance. Our ability to reliably measure the parameters of a binary SMBH depends crucially 
on the signal model onto which the data is projected. If the signal model is incorrect or even inaccurate, then 
regardless of how good the instrument might be, the measurement is prone to systematic errors which will make 
precision measurements meaningless (sec, e.g., Ref. [45j for a detailed discussion on parameter extraction errors due 
to inaccurate template waveforms) . From the point of view of incurring the minimum possible errors in the estimation 
of parameters it is desirable to use the best known signal model. We have provided critical insight into the role of 
higher harmonics in the fast gravitational-wave phasing, higher FN corrections in the amplitudes and their interplay 
with the slow modulations induced by LISA's motion (see Appendix B). 

By using a waveform whose phase is correct to order (3.5PN) and with the amplitudes correct to order (2.5PN), 
we have found that LISA's angular resolution improves typically by a factor of 10, even with a very conservative choice 
for the detector's lower cutoff frequency (i.e., 10~^ Hz). This means that LISA will be able to uniquely identify the 
galaxy cluster in which the merger event took place and thereby facilitate optical identification and the measurement 
of the redshift of the source. Together with the fact that binary SMBH sources are standard sirens, this means that 
LISA will be able to measure the cosmological parameters by circumventing the lower rungs of the cosmic distance 
ladder. 
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APPENDIX A: IMPLEMENTATION OF THE PARAMETER ESTIMATION WITH FWF AND ISSUES 

RELATED TO NUMERICAL ACCURACY 

The computation of the covariance matrix for LISA using FWF starts with a simpler code (in Mathematica) which 
reproduces the results of Ref. [H, [13, [2^ for the ground-based detectors. The modulations due to LISA's orbital 
motion are then accounted for to obtain the waveform in LISA's frame. This waveform is differentiated analytically 
and the derivatives are stored in an array and used in the computation of the various elements of the Fisher matrix. 
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The integrations that are needed in computing different elements of the Fisher information matrix are performed using 
numerical integration routines of Mathematica. The program thus evaluates the Fisher matrix for the given input 
values of the component masses of the binary, the four angles describing the source in the solar system barycenter 
frame and the PN order of the amplitude for a fixed distance of 3 Gpc (phasing is always 3.5PN). The input masses 
are always physical masses consistent with the cosmological model we described earlier. 

We have checked that our code reproduces the numbers in Ref. [IHl (with the corresponding restricted waveform 
which uses a 1.5PN phasing including the spin-orbit parameters, and the corresponding noise PSD) with less than 
1% difference. Since we use the noise PSD of Ref. Q, we have checked that we obtain the numbers for the pattern- 
averaged case of to less than 1% difFercncc using their waveform parametrization which includes all nonspinning 
and spinning binary parameters up to 2PN but ignores higher PN order terms in the phasing. Lastly, we recover the 
results of [231, which uses a 3.5PN phasing for the parameter estimation and recover the results again to less than 1% 
difference. 

The numbers quoted in this paper were checked with three independent codes, all in Mathematica, but using 
different numerical integration routines (such as NIntegrate and Listlntegrate) and matrix inversion routines 
(such as Inverse and another method based on singular value decomposition (SVD)), and all three codes agree to 
less than 5%. 

A large dimension of the parameter space often leads to ill-conditioned Fisher matrices, i.e., the magnitude of 
the ratio of the smallest and the largest of its singular values (the inverse of the condition number) approaches the 
machine's floating point precision (10^^^ in our case). For the (6.45 x 10'', 1.29 x 10^)Mq system, only the RWF 
Fisher matrix is ill-conditioned, for the third angle. The (1O^,1O^)M0 system, for both RWF and FWF, is free of 
ill-conditioned Fisher matrices. More specifically, the inverse of the condition number is always 10^ — 10'^ times larger 
than 10~^^. The Mathematica inversion routine as well as singular value decomposition (SVD) are used to obtain 
the covariance matrix and both methods give the same results. For these systems, we also note that the numerical 
eigenvalues (computed using the Mathematica function Eigenvalues) coincide with the numerical singular values 
(obtained through SingularValueDecomposition) at the standard floating point precision of 10~^^. On the other 
hand, the Fisher matrices for the (10®, 10^)Mo system show ill-conditioned behaviour at this precision. However, this 
may not mean that the results obtained are unreliable for reasons next outlined. The SVD of the Fisher matrices 
shows that one or more of the singular values approach zero. Indeed, a machine precision calculation of the singular 
values yields zero for at least one singular value of the Fisher matrices for this system. In such cases, we replace the 
singular values by the numerical eigenvalues (whose machine precision calculation does not yield zero) to obtain the 
inverse using SVD. However, Mathematica has the option of computing at higher precision by padding unknown 
digits beyond the known ones. Using this facility and repeating the inversion procedure using SVD but with precision 
higher than 10~^® does not alter the results and shows that the singular values are equal to the eigenvalues, as 
expected, with none of them being zero. The condition number is also seen to be 10 times larger than the precision 
used. In the case of one or more almost zero singular values (at the standard floating point precision), one can also 
obtain a pseudoinverse [i^ that is closest to the "real" inverse, in the least-squares sense, and end up with different 
results. However, for the (10®, 10^)Mq system, the errors given by the pseudoinverse do not seem physically correct 
and so we refrain from using it in this work. A detailed analysis, by perturbing the Fisher matrices and observing 
how the inverses behave will be a stronger test of the reliability of our results [43| . 

A recent Markov Chain Monte Carlo (MCMC) analysis of Ref. [i^ revealed an interesting point relevant to our 
analysis. It compared the estimates from the Fisher information matrix with those obtained from MCMC and found 
that there is excellent agreement between the two methods in the case of the extrinsic parameters (what we refer to 
as Psiow), although in the case of the intrinsic parameters the two methods are not in good agreement. This is good 
news since our main concern here is the angular resolution and the luminosity distance. Though their analysis used 
the restricted waveform, we believe that similar trends would exist in the case of FWF too since the dimensionality 
of the parameter space is the same as that of RWF. 

APPENDIX B: PARAMETER ESTIMATION WITH THE FWF AND VARIATION WITH PN ORDER IN 

AMPLITUDE 

In this Appendix we will attempt to provide some qualitative understanding of the trends observed in the parameter 
estimation accuracy with increasing PN order in amplitude. Recall that the covariance matrix goes as 1/SNR^ and 
therefore errors go down as 1/SNR. However, notice from Table |T] that even though the FWF SNR is less than the 
RWF SNR all errors decrease as we go from RWF to FWF. This provides the clue that the SNR cannot account for 
the improved performance of FWF: In going from the RWF to 0.5PN, the ratios of the two SNRs are nowhere close to 
the inverse ratios of the errors in various parameters. Further, as we go from 0.5PN to 2.5PN, although the products 
of the SNR with errors do not fluctuate much, they are noticeably different implying that the SNR is not the main 
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90 


65 


32 


2400 
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0.5 
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32 


2.6 


2.3 


6.3 


9.9 


0.5 (3rd) 


2004 


110 


110 


32 


2.6 


2.3 


6.3 


10 


0.5 (1st) 


1825 


370 


310 


81 


21 


10 


360 


960 


0.5 (2nd) 


1648 


410 


340 


90 


12 


10 


690 


1800 


0.5 (2nd with span of 3rd) 


1757 


190 


140 


31 


5.0 


2.2 


17 


44 


0.5 (7th) 


1943 


11 


22 


8.5 


1.9 


2.1 


3.6 


4.2 



TABLE III: Parameter estimation at 0.5PN order with mock waveforms corresponding to different choices. The first (second) 
row corresponds to the restricted (0.5PN) results. The third (fourth) row corresponds to the 0.5PN waveform with only the 
3rd (1st) harmonic and suppressing the 1st (3rd) harmonic. The fifth row corresponds to a mock waveform suppressing the 
first harmonic and replacing the third harmonic by the second. This implies retaining its polarisation amplitude and phase 
but changing f(//3) — * i(//2) in the Fourier and the Doppler phases and decreasing the frequency span from 3(Flso — ^"in) 
to 2(_Flso — -Fin). The sixth row corresponds to the previous case but with an enhanced span 3(-Flso — ^In). The last row 
corresponds to a mock waveform where the seventh harmonic is moved from 2.5PN to 0.5PN and treated with its associated 
normal frequency span 7(-fLso — ^In) with the 1st and 3rd harmonics suppressed. In this and other tables in this appendix the 
system and angles are the same as in Table U 

reason why the parameter estimation improves. 

Next, we note that at 0.5PN both the first and the third harmonic possess the same structure (/ dependence) and 
almost the same variety (dependence on the inclination angle i) (compare Eq. (5.7b) and (5.8b) of [2J|). However, 
their effects on PE, when examined independently, are vastly different as seen in Table Hill The 0.5PN PE is hardly 
distinguishable from the PE with the first harmonic suppressed and thus the improvement at 0.5PN is solely due to 
the third harmonic. This is further confirmed by performing the PE with a "mock" FWF where the third harmonic 
is replaced by the second harmonic: as seen from Table [TTTl the PE worsens in this case. Thus, the presence of higher 
harmonics brings about improved PE. It must, however, be pointed out that the improvement is obtained not merely 
due to the presence of the higher harmonic but crucially due to the increased span of the fcth harmonic in the frequency 
domain (A;(-FLSO~-Pin))- This is first checked by working with a mock 0.5PN waveform containing the second harmonic 
but using for its span the increased span of the third. It is confirmed by an alternative mock waveform, suppressing 
the first and third harmonics and including instead the seventh harmonic at 0.5PN. The dramatic improvement in 
PE compared to the regular 0.5PN waveform, but with almost the same SNR, shows that at the 0.5PN level, where 
complications due to PN corrections to the harmonics and higher order expansion coefficients of F~^l'^ arc not present, 
the SNR is not the determining factor. 

To understand the variation in PE from 0.5PN to 2.5PN we need to disentangle many effects. The FWF introduces 
more structure (/ dependences) at different PN orders. Different harmonics appear at different PN orders and in our 
model are associated with their respective spans. At 2.5PN, higher order corrections appear in every harmonic except 
the sixth and the seventh harmonics. The frequency sweep F~^l'^ in the Fourier transform of the binary signal brings 
in further PN corrections associated with it (the coefficients S'(m/2) in Eq. (|2.3p ). Finally, the coefficients involve angles 
which are fixed in the barycentric frame but suffer modulations (Doppler and orientational) due to LISA's motion. 
To this end, we reexamined the PE with a series of mock waveforms including where possible only one of the above 
aspects. For instance, in Tabic [TVl we consider only the second harmonic and its associated PN corrections. From the 
results it is clear that although PN corrections to the 2nd harmonic bring in additional terms, they do not improve 
PE; rather, in comparison with the RWF, the PE worsens. This is probably due to the fact that PN corrections with 
the same frequency dependence as before (or powers thereof) increase the covariances among the various parameters 
and thus worsen PE. 

In Table IVl on the other hand, the mock waveform includes all the higher harmonics but excludes their PN amplitude 
corrections and PN corrections arising from the frequency sweep. Clearly, the PE improves monotonically as we go 
from one order to the next in all the parameters quoted and so docs the SNR. Additionally, an examination of the 
results in Tables [I] and |V] reveals that while the higher harmonics by themselves tend to generically improve PE, 
the PN corrections to the harmonics arising from the higher PN amplitudes and higher PN frequency sweep tend to 
degrade PE. However, their effect is less dramatic than one may naively expect. This is because higher harmonics 
necessarily appear at higher PN orders and, as is evident from Eq. (|2.3p . they come with higher powers of the PN 
expansion parameter: (27rAf fk)"'^^ for the fcth harmonic appearing at the n/2th PN order. The upper cutoff of the fcth 
harmonic in the frequency domain is fci^LSO a-t which (27rAf/)^/'^ reaches its maximum of 6~^^^. The 7th harmonic. 
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1353 


500 


410 


120 


130 


48 


2400 


7200 



TABLE IV: Variation of PE with PN orders in amplitude for a mock FWF retaining only the 2nd harmonic and its higher 
order PN corrections. The PN corrections to a given harmonic at PN order n/2 add terms of type {2nM f /k)"^^ . While they 
bring in new structure in the waveform they do not help improve PE; instead they enhance the covariances among different 
parameters and thereby worsen PE relative to RWF. 
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TABLE V: Variation of PE with PN orders for a mock FWF retaining all the harmonics and neglecting both higher order 
PN corrections to them and also PN corrections to i^)^^^'^ arising from the frequency sweep. Thus most higher order PN 
corrections to the harmonics are neglected. Higher harmonics generally improve parameter estimation. PN corrections to 
harmonics tend to degrade PE. 

for example, will be scaled by a factor which is always less than 6^^^^ and consequently will contribute less power. 
The seventh harmonic in the FWF at its regular PN order hardly leads to any improvement in PE. However, as 
remarked earlier in Table IIIIl had it appeared at 0.5PN its impact would have been substantial. We have verified 
these trends for another angle and for a lower mass system of (10^, 10^)A/q; for the sake of brevity we leave out the 
details. 



(10^ 10^)Afo; Correlation- Coefficient Matrix; OPN (A3) 




\nM 


S 


tc 


0c 


Ini^L 


MS 


0s 


ML 




InM 


1 


-0.99 


0.93 


-0.19 


-0.17 


0.14 


-0.13 


-0.15 


0.031 


S 




1 


-0.92 


0.098 


0.084 


-0.06 


0.065 


0.074 


-0.0088 


tc 






1 


-0.43 


-0.45 


0.33 


-0.43 


-0.46 


-0.022 


0c 








1 


0.98 


-0.95 


0.70 


0.92 


-0.39 


lnL»L 










1 


-0.88 


0.81 


0.98 


-0.22 


MS 












1 


-0.45 


-0.77 


0.65 


0s 














1 


0.92 


0.38 


ML 
















1 


-0.016 


0L 


















1 



TABLE VL Correlation coefficients computed using RWF for a binary at 2 ~ 0.55 comprising (10*^, 1O'^)M0 SMBH. The Table 
shows how the 'fast' and 'slow' variables behave. There is a high correlation among parameters in the same subclass, but only 
a weak correlation between members of different subclasses. Entries that are vacant can be found by symmetry. 

Our last comment relates to the decoupling between the fast variables pfast and the slow variables Psiow in the case 
of the FWF. We have checked that in the case of the FWF the correlation coefficients are very small (of the order 
of 0.1-0.3) between these two subsets of parameters (see Table IVTl and Table |VlT|. There exists, however, a high 
correlation among different parameters of the same type as in the case of RWF (sec, e.g., Ref. @). Specifically, as in 
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(10^ 10'^)A/o; Correlation-Coefficient Matrix; 2.5PN (A3) 
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TABLE VII: Correlation coefficients computed using FWF for a binary at 2 ~ 0.55, comprising (10*^, 10^)Mq SMBH. The 
Table clearly shows a more pronounced decoupling between the 'fast' and 'slow' variables at 2.5PN. There is a high correlation 
among elements of the same subclass but a lower correlation between the fast and slow subclasses. 

the RWF case , the distance estimation could be improved if the source is better localized in the sky. 

It is worth recalling that in the parameter estimation problem we are not only limited by statistical errors due to 
noise but also by theoretical or s yste matic errors arising at any PN order due to the limited accuracy of the waveforms. 
Recently Cutler and Vallisneri [4^ have looked into this issue more critically. In the context of the present work, 
we would like to stress that despite the fact that the largest improvement in PE arises from the third harmonic at 
0.5PN, the need to limit systematic errors mandates the use of the best available waveform, i.e., the one at 2.5PN in 
amplitude and 3.5PN in phase. 



